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PACS 05 . 45 . -a - Nonlinear dynamics and chaos 

Abstract. - We generalize the concept of convective (or velocity-dependent) Lyapunov exponent 
A(ii) to an entire spectrum A(v,n). Our results are supported by the consistency between the 
outcome of the chronotopic approach [S. Lepri et al. J. Stat. Phys., 82 5/6 (1996) 1429] and 
a more direct method. There exists a critical integrated density n — n c , beyond which the 
convective exponent exhibits a discontinuous dependence on the velocity, which originates from 
the appearance of multiple branches. This phenomenon can be traced back to a change of concavity 
of the so-called temporal Lyapunov spectrum for n > n c , which is therefore a dynamical invariant. 
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Introduction. — The linear-stability analysis of an 
N dimensional chaotic dynamics typically amounts to 
computing N Lyapunov exponents (LEs) {A;} with I = 
1, . . . , N (typically ordered from the largest to the most 
negative one). In a spatially extended system, the Lya- 
punov spectrum depends only on the integrated density 
n = l/L (A(n)) if the system size L is large enough pQ. 
Since this tool does not provide any information on the 
spatial propagation of perturbations, a new indicator was 
introduced in Ref. [2J, the convective (comoving) Lya- 
punov exponent (CLE) A(u), which quantifies the max- 
imal growth rate of an initially localized perturbation, in 
a frame that moves with a velocity v. Later, it was rec- 
ognized that the CLE can be derived from a general the- 
ory, based on the so-called chronotopic approach [3] , which 
deals with the wider class of perturbations with a spatially 
exponential profile exp(/ii), where i denotes a discrete spa- 
tial variable and /i is a free parameter. As a result, one 
can, e.g., define the generalized temporal Lyapunov A(n, fi) 
spectrum. The CLE A(v) can be thereof computed as a 
Legendre transform of the maximal temporal Lyapunov 
exponent, i.e. A(0,/x) [4]. 

In systems with left-right spatial symmetry, A(v) is 
symmetric around v = 0, where it attains its maximum 
value which coincides with the standard maximum LE. 
The largest velocity v = v c such that A(v) > is the max- 



imal propagation velocity of infinitesimal perturbations. 
In convectively unstable systems A(0) < and the (posi- 
tive) maximum convective exponent is attained for some 
nonzero velocity. In fact, a somehow similar approach was 
developed by Huerre and Monkewitz [5] to characterize ab- 
solute and convective instabilities in open flows and more 
recently extended by Sandstcde and Schecl to deal with 
generic boundary conditions [5J. 

In this Letter we go beyond the computation of the 
growth rate of the local amplitude of a perturbation, turn- 
ing our attention to the evolution of volumes of generic 
dimension. The chronotopic approach offers a straightfor- 
ward way to determine an entire spectrum of convective 
LEs, by extending the notion of Legendre transform from 
the maximum exponent (i.e. n = 0) to generic values of 
the integrated density n. However, it would be desirable 
to give a more direct definition as well. In principle, the 
most appropriate starting point for a direct definition of a 
spectrum of CLE is the approach developed in [7], where 
an ensemble of linearly independent initial conditions (lo- 
calized within a window of size L) was freely let evolve 
to thereby determine volume expansion rates within the 
very same window. As a result, it was noticed that mean- 
ingful generalized Lyapunov exponents could be defined 
by simultaneously letting the time T and the size L tend 
to infinity, with constant ratio g = T/L. The standard 
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Lyapunov spectrum is recovered for g — Y 0, while in the 
opposite limit g — > oo, all exponents coincide with the 
maximum. The dependence on g expresses the fact that 
the growth along different directions is affected in a dif- 
ferent manner by the local expansion and by diffusion. 
In principle, this approach can be implemented for mov- 
ing windows, too; however, the unavoidable dependence of 
the additional parameter g, reduces the appeal of such a 
direct method. Accordingly, we have preferred to comple- 
ment the moving-window approach with some boundary 
conditions to get rid of the difficulty of dealing with an 
open system. In practice, this is the setup proposed in 
Ref. 0. 

Altogether, in this Letter we compare this latter method 
with the chronotopic approach, finding that they are mu- 
tually consistent. Moreover, the chronotopic approach is 
by and large the most accurate and this has allowed dis- 
covering serious numerical difficulties that easily lead to 
artifacts in the direct computation of the convective Lya- 
punov spectra (even in the simple case of coupled maps, 
herein investigated). Moreover, the chronotopic approach 
has revealed that the lower part of the convective spectrum 
is charactrized by the existence of multiple solutions. The 
phenomenon appears for n > n c and is associated to a 
change of concavity in the temporal Lyapunov spectrum. 
The critical density n c is, by construction, a dynamically 
invariant dimension density (analogous to the dimension 
density of the stable manifold or to the Kaplan- Yorkc 
dimension density), but the physical meaning of n c still 
needs to be clarified. 

Altogether the consistency between the two methods 
confirms the conjecture that the chronotopic approach 
"encodes" all stability properties of one-dimensional 
spatio-temporal systems, that are eventually contained in 
the corresponding entropy potential j^j. 

Theory. — First, we introduce the proper formalism 
with reference to a standard model of coupled maps, 



Vl+i = f 



(l-e)yi + -{y\- x +yi +l ) 



(1) 



where i = 1, . . . , N and t are the spatial and temporal (in- 
teger) indices, respectively, while e £ (0, 1) represents the 
diffusive coupling and f(y) is a map of the unit interval 
onto itself. We start by investigating the evolution of an 
infinitesimal perturbation 8y\ } initially localized in a finite 
region of length M + 1 < < N centered around the origin 
(Syi = CQ(M/2 - i)Q(i + A//2), where 9 is the Heavi- 
side function and the £ l 's are iid random variables). The 
standard convective Lyapunov exponent is defined as 
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where the "0" in parentheses signifies that we refer to the 
maximum exponent. The formula can be easily extended 
to generic dynamical systems, by replacing the absolute 
value \5yl\ with any norm quantifying the amplitude of 



the perturbation on the site i at time t. The question 
addressed in this Letter is the extension of this definition 
to quantify not only the growth rate of the amplitude but 
also the expansion/contractions along the additional di- 
rections. In order to address this issue, it is necessary 
to explore the spatial structure of the perturbations. In 
the standard approach, it is sufficient to first let the per- 
turbation evolve freely in the space and, afterwards, to 
focus the attention on specific world lines (i = vt). Will- 
ing to go beyond the maximal convective LE, one should 
consider a set of localized but linearly independent initial 
conditions. Accordingly, a theoretical and practical diffi- 
culty emerges: that of defining a proper orthogonalization 
strategy to avoid the convergence towards the same di- 
rection. In the case of the standard Lyapunov spectrum, 
this problem is solved by resorting to the Gram-Schmidt 
technique [p_j. However, here it is not obvious how to do 
that. The natural idea of referring to the space covered 
by the perturbations at a given time leads to the unavoid- 
able conclusion that one would measure only the tip of 
the spectrum, since we would have to rescale by a growing 
coefficient. Therefore, one is led to restrict the orthogo- 
nalization to a moving window of fixed size L < M + 1, so 
that the computation of convective Lyapunov spectra for 
different velocities would require performing separate sim- 
ulations. This idea, which was already suggested in [8], is 
not the optimal solution, as the initial perturbations are 
not let to freely evolve, but artificial modifications have 
to be introduced on the boundaries to confine the evo- 
lution within the prescribed window. Nevertheless, since 
it turns out that boundary conditions do not matter for 
large enough L (see below) , we can at least claim that the 
approach is meaningful. 

Let us now be more specific and explain how the LE can 
be determined in a frame that moves with a generic veloc- 
ity v in a discrete spatio-temporal lattice. We introduce 
two types of "moves" in tangent space: the first corre- 
sponds to a static window (0); the second to a right-shift 
by a single site (1, without loss of generality we discuss 
only windows moving to the right). The corresponding 
rules are, 



Syt +1 = ml 



-Syr 1 + (I - e)8yl + -Syf +1 



(3) 



where 1 < k < L, i denotes the absolute position in the 
lattice and m\ = /'[(I -e)y\ +e{y\~ 1 +y\ +1 )/2]. 0- and 1- 
iterations correspond to j = k and j = k + 1, respectively. 
The restriction of the rule to a finite interval requires ex- 
tra assumptions for 5y^ +1 and Sy® (Sy^ +2 ), to close the 
model for a (1) move. We have typically worked by 
assuming all of them to be equal to zero, but we have ini- 
tially verified that the same results are obtained also for 
different choices (e.g., Sy® = 8y\). As a last detail, it is 
necessary to specify the sequence of 0s and Is that is used 
to study the velocity v. Given that v is obviously equal to 
the fraction of Is, we have typically selected (simple) ra- 
tional numbers, choosing the most uniform sequence with 
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a fixed density (e.g., for v = 2/5, 01010010101001 . . .). 

As a result, for any given velocity v, we can determine the 
whole spectrum A(v,n), where n = l/L and I means we 
refer to the Ith largest exponent. 

An alternative approach consists in considering the 
evolution of a perturbation with an exponential profile, 
namely 8y\ = $Je'". This requires studying the evolution 
equation 



t+i 



(1 



(4) 



where i = 1, . . . , N and periodic boundary conditions for 
$J can be safely assumed. By solving the model (Q| for 
different values of fx, one obtains the temporal spectrum 
X(n,n). 

We claim that the convective Lyapunov spectrum can 
be obtained by a Legendre transform of A(n, /x), 



A(i/, n) = A(n, fi) — fiv 



(5) 



where v = dX/dfj, (and the set of transformations is com- 
pleted by the symmetric expression fi = dA/dv). Alto- 
gether Eq. ([5]) generalizes to n > the transformation 
that has been proved for n = [1]. Its original justifi- 
cation is based on the observation that the profile of an 
initially localized perturbation is, at time t, locally expo- 
nential with a decay rate \i that depends on the position 
i = vt. For 7i > 0, it is not clear which profile one should 
refer to. For this reason the above definition is essentially 
a formal one. 

Numerical Results. — We have tested both ap- 
proaches, by studying a chain of coupled logistic maps, 
i.e. with f(y) = ry(l — y) (r = 4), y G [0; 1] and e = 1/3 
(these are the same parameter values adopted in Ref. [5])- 

In Fig. Q]wc plot the temporal spectra for N — 100 and 
some values of n (namely, n = 0, 0.5, 0.75, 0.9 and 1). 
We have also verified that N is large enough to ensure 
the thermodynamic limit. The uppermost curve corre- 
sponds to the maximum exponent that is known to grow 
monotonously for increasing (in absolute value) \i. On the 
other hand, the lowermost curve, which corresponds to 
the minimum LE not only is non monotonous, but even 
exhibits a singular behavior for ?» 1.31. 

Convective exponents vs velocity. In Fig. [2^, we plot 
the convective spectra obtained by Legendre transform- 
ing two of the curves reported in Fig. [TJ The upper solid 
curve is the standard convective spectrum. It starts from 
a maximum value for v = that corresponds to the usual 
maximum LE and crosses the zero axis at v c = 0.510(5) 
which indicates the maximal propagation velocity of per- 
turbations. The lower solid curve corresponds to n = 0.5 
and the very fact it is nearly zero for v = indicates that 
approximately half of the standard LEs arc positive. 

In order to check the equivalence between the direct 
and the chronotopic approach, we have computed the con- 
vective spectra also by iterating localized perturbations, 




Fig. 1: (Color Online) Temporal Lyapunov exponents A(n, /i) 
versus /j for a chain of iV = 100 coupled logistic maps with 
r — 4 and e = 1/3. From top to bottom the curves refer to 
ra = 0, 0.5, 0.75, 0.9 and 1. 
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Fig. 2: (Color Online) Convective exponents A(v, n) vs. the 
velocity of the comoving reference frame for Logistic (a) and 
Bernoulli (b) coupled map lattices, (a) The convective expo- 
nents have been obtained for n = (upper curve) and n = 0.5 
(lower curve). The solid lines have been obtained be Legendre 
transforming (see Eq. [5| the curves reported in the previous 
figure. The symbols correspond to direct estimates of the con- 
vective exponents in double precision for n = 0.5 and L — 100: 
(orange) asterisks refer to 7 = 0, while (red) squares to 7 = 1.0. 
(b) The reported convective spectrum has been analytically de- 
termined by performing the Legendre transform of Eq. © for 
n = 0.75 and N = 100. 
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implementing the method described in the previous sec- 
tion. The (orange) asterisks in Fig. correspond to the 
outcome of simulations performed with L = 100. The 
agreement is rather poor and does even get worse upon 
increasing the window size L (not reported). However, 
we discovered that the problem is not a conceptual, but a 
numerical one. The reason is that the Lyapunov vectors 
in the moving windows are exponentially localized around 
the left border of the window. Accordingly, many of their 
components are so small that double-precision accuracy 
(i.e. with 14-15 digits) is not sufficient, especially for large 
values of L. As a matter of fact we have verified that by 
employing extended precision (i.e. with approximately 30 
digits) the agreement increases. However, an even better 
agreement can be obtained already with double precision, 
by using weighted Euclidean. More precisely, given any 
two vectors u = {u 1 }, v = {V}, we define the scalar prod- 
uct as 

u • v = ^2 UiVie 11 (6) 

i 

where 7 is a free parameter. Formally speaking, it is well 
known that Lyaphnov exponents are independent of the 
norm that is chosen, i.e., in this case, of 7. Nonethe- 
less, the numerical accuracy may significantly depend on 
7. In fact, by comparing double with extended precision 
for different values of 7 we concluded that 7 = 1 is a 
nearly optimal choice Q. The data reported in Fig. [2^i (red 
squares) indeed confirm the increased agreement with the 
chronotopic results. However, a proper selection of 7 does 
not solve completely the problem: for large velocities and 
larger values of L, numerical accuracy remains a serious 
issue. 

Anyway, the relevant message that comes from the nu- 
merical analysis is that the definition of convective Lya- 
punov spectra through the chronotopic approach is not 
just formal but provides the correct answer and, more- 
over, the method is far more reliable than the direct one, 
since one does not face accuracy problems. Incidentally, 
it is the comparison between the two approaches that has 
allowed discovering the serious problems of numerical ac- 
curacy that affect the direct approach and thereby the sim- 
ulations reported in Ref. |S|. Last but not least, for large 
values of the integrated density n, the non-monotonicity 
of X(n, /i) leads to the emergence of extra branches in the 
convective spectrum, as shown in Fig. [2b (which refers to 
Bernoulli maps f(x) — 2x mod(l) with the same coupling 
constant as before). The origin of the three branches at 
small velocities will be addressed in the following section, 
while the coexistence of branches at large velocities is to- 
tally irrelevant. In fact, it arises in an unphysical region, 
since v > 1 cannot be obtained in a lattice with nearest 
neighbor coupling. We interpret the very existence of the 
two branches as an instance of a phase velocity that can be 
faster than "light velocity" without causing any paradox. 

1 At least for not too large velocities, otherwise the choice of the 
norm does not help enough 
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Fig. 3: (Color Online) Spectrum of the convective exponents 
A(v,n). (a) Coupled logistic maps for v — 1/5. The solid lines 
correspond to direct measurements. From top to bottom: L — 
100 and 7 = (red dashed line); L = 200 and 7 = 2 (blue dot- 
dashed line); L = 100 and 7 = 2 (green solid line). Full circles 
correspond to the Legendre transform for N = 100. Finally, 
the vertical dot-dashed (magenta) line indicates the n c -value. 
(b) Bernoulli maps for v — 1/3, estimated analytically via the 
Legendre transform (red circles) and directly by diagonalizing 
a constant 50 x 50 matrix (solid line) with extended precision. 

Convective exponents vs density. We now compare the 
convective spectra obtained for a given velocity, as it helps 
clarifying the behavior for large integrated densities. An 
example is reported in Fig. [3^ for v = 1/5. Above a cer- 
tain density n c the spectrum obtained via the Legendre 
transform displays three different branches: the lower one 
is associated to positive /^-values, while the other two to 
negative \i. Once again the standard direct estimates suf- 
fer problems of numerical accuracy (see the upper solid 
curve in Fig. ^l). The improved simulations (for 7 = 2) 
seem to converge towards the lower branch, but there is 
still some discrepancy and it is difficult to decide whether 
this is due to finite-size corrections. 

In order to obtain a more convincing evidence, we now 
consider a model of Bernoulli maps. Since the multipli- 
ers are constant in space and time, we expect smaller fi- 
nite size corrections and, moreover, the estimation of the 
convective spectra is simpler as it reduces to the diago- 
nalization of a matrix (for rational velocities). In fact, 
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Fig. 4: (Color Online) Isolines in the plane (A, fi) corresponding 
to a constant velocity v for coupled logistic maps. The vertical 
dotted line refer to v = 0, the other dashed lines from left to 
right for /i > (resp. right to left for fi < 0) correspond to 
v = 0.05, 0.20, 0.50 and 0.95 with N = 100. 

in this model, the tangent matrices (see Eq. depend 
only whether the move is of type or 1. Therefore, given 
a rational velocity, characterized by a periodic sequence of 
Os or Is (with period P), it is sufficient to multiply P of 
such matrices and thereby determine the eigenvalues. The 
results for v = 1/3 = (001001001 . . .) are reported in Fig. 

where one can appreciate that the agreement between 
the lower branch and the direct method is already quite 
impressive for a window of size L = 50. 

In order to shed some light on the origin of these 
branches, we have plotted the paths in (A,^i)-plane that 
correspond to spectra with different velocities (sec Fig. @| . 
They are the isolines where the slope dX/dfi = v stays 
constant. The vertical line at fx = is the path for the 
standard Lyapunov spectrum. However, for finite veloc- 
ities the line breaks into two components that lie in the 
positive and negative fi half-planes, respectively: the for- 
mer one corresponds to the lower branch in Fig. [3J while 
the latter one gives rise to the two upper branches. Al- 
together it is reasonable to conjecture that the reason for 
discarding the upper branches is that they correspond to 
negative fi values, i.e. to profiles that are larger on the 
right side. In fact, this is impossible unless one assumes 
the presence of some source of "noise" on the right of the 
moving window. 

A moment's reflection suggests that the mathematical 
origin of two branches in Fig. [4] is the change of concav- 
ity of A(n,/x) (as a function of n) upon increasing n: it 
suddenly enforces points with positive derivative to jump 
from the right to the left side of fj, = 0. Once again the 
Bernoulli maps allow for an accurate investigation. 

The expression for the temporal exponent is [3], 

A(n,/i) = logo + - log |(1- ef+ (7) 
2e(l — e) cosh \x cos 7m + e 2 (cosh 2 fi — sin 2 7m) | 
By expanding this expression for small values of /x, it is 
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Fig. 5: (Color Online) Spectrum of the convective exponents 
A(v, n) for coupled Bernoulli maps for v = (circles) and v = 
0.001 (squares), estimated by Legendre transforming of Eq. (0). 
The parameters are a — 2 and e = 1/3. 

easy to verify that the quadratic term in \x changes sign 
for 

_ acos[e/(l - e)} 
n = n c = l , (8) 

7T 

provided that e < 1/2. As a matter of fact we have esti- 
mated the spectrum for the convective exponents for v = 
and v = 0.001. In Fig. [5j one can appreciate that for 
n < n c they are nearly identical, while for n > n c , they 
separate out by a finite amount as a result of the selection 
of the lower branch. Now, we can return to A(v,n) as a 
function of v for fixed n (illustrated in Fig. [2] for n = 0.5) 
and conclude that A(u, n) must have a discontinuity in 
v = 0, if n > n c . 

Discussion and conclusions. In this Letter, we 
have shown that the notion of convective Lyapunov ex- 
ponent can be extended from the maximum to an entire 
spectrum. The comparison between two different methods 
reveals a rather complex scenario with several issues that 
need be further clarified. An example is the dimension 
density n c which separates a standard behavior from the 
appearance of multiple branches. Its very existence is con- 
nected to a change of concavity in the temporal Lyapunov 
spectra. We have verified that in coupled maps it gener- 
ally exists for not-too-largc coupling (e < 1/2 in Bernoulli 
maps) and preliminary simulations confirm the existence 
of such a critical density also in Stuart-Landau oscillators 
[ID] . Nevertheless, its physical meaning is rather unclear. 
We can only claim that the dimension n c is a dynamical 
invariant (as it follows from a general property of the Lya- 
punov spectrum) like the density of unstable directions 
(fraction of positive Lyapunov exponents), the Kaplan- 
Yorke dimension and the dimension of physical modes |11] . 
Another question concerns the upper branches that we 
have dismissed as irrelevant, but could play some role in 
specific physical contexts. 

Finally, the direct definition is still partially unsatisfac- 
tory, as it involves the addition of in-principle-unnecessary 
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boundary conditions. An open system approach such as 
that developed in [7] would be much more appealing, but it 
requires incorporating the additional parameter g = T / L 
in the current theory. 

An alternative idea could be that of referring to covari- 
ant Lyapunov vectors [H] : if one could indeed "build" the 
initial perturbation by using only such vectors, one would 
be automatically assured that no more-unstable degrees 
of freedom are going to be excited. The problem is that 
one first needs to evolve back and forwards the chain on 
sufficiently long time scales to allow the perturbation prop- 
agate over sufficiently large distances to measure asymp- 
totic quantities. 
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